function [ P, output ] = FrankWolfe( D, lambda, para )
% this version need to know testing points

if(~isfield(para, 'pwIter'))
	para.pwIter = 15;
end

tol = para.tol;
maxIter = para.maxIter;

[row, col, data] = find(D);

P = zeros(length(para.row), 1);
X = zeros(length(row), 1);
spa = sparse(row, col, data, size(D,1), size(D,2));

obj = zeros(maxIter, 1);
RMSE = zeros(maxIter, 1);
v = randn(size(D,2), 1);
for i = 1:maxIter
    s = X - data;
    obj(i) = (1/2)*sum(s.^2);
    
    % approximate SVD
    spa = setSval(spa, s, length(s));
    [u, pwIter] = powerMethod(spa, v, para.pwIter, 1e-3);
    v = -(u'*spa);
    v = v'/sqrt(sum(v.^2));
    
    % update gradient
    gamma = 2/(i + 1);    
    X = (1 - gamma)*X +(gamma*lambda)*partXY(u', v', row, col, length(row))';
    % update prediction
    P = (1 - gamma)*P +(gamma*lambda)*partXY(u', v', para.row, para.col, length(para.row))';
    
    if(i < 6)
        delta = inf;
    else
        delta = abs(obj(i-5:i-1) - obj(i-4:i));
        delta = mean(delta);
    end
    if(delta < tol)
        break;
    end
    if(mod(i + 1, 100) == 1)
        fprintf('iter:%d obj: %.3d dif:%.2d, pwiter %d \n', ...
        i, obj(i), delta, pwIter);
    end
    if(isfield(para, 'test'))
        RMSE(i) = sqrt(sum((P - para.test.data).^2)/length(P));
    end
end
output.obj = obj(1:i);
output.RMSE = RMSE(1:i - 1);

end

